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Abstract 

We consider heat transport in one-dimensional harmonic chains 
with isotopic disorder, focussing our attention mainly on how disorder 
correlations affect heat conduction. Our approach reveals that long- 
range correlations can change the number of low-frequency extended 
states. As a result, with a proper choice of correlations one can control 
how the conductivity k scales with the chain length N. We present 
a detailed analysis of the role of specihc long-range correlations for 
which a size-independent conductivity is exactly recovered in the case 
of hxed boundary conditions. As for free boundary conditions, we 
show that disorder correlations can lead to a conductivity scaling as 
K ~ , with the scaling exponent e being arbitrarily small (although 

not strictly zero), so that normal conduction is almost recovered even 
in this case. 
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1 Introduction 


The derivation of the transport properties of one-dimensional (ID) models 
from the underlying microscopic dynamics is a long-standing problem in non- 
equilibrium statistical mechanics. An especially relevant (and vexed) ques¬ 
tion is: What conditions are required for heat conduction to obey Fourier’s 
law (see, e.g., [1])? Among the various ID models, harmonic chains are the 
most convenient for analytical treatment and have been intensively studied 
over many decades. The protracted research efforts have led to the conclusion 
that Fourier’s law does not hold in ID harmonic chains with isotopic uncor¬ 
related disorder [i, 2, 3], unless the baths are endowed with special spectral 
properties [4]. The interest in this topic was signihcantly heightened after 
experiments provided evidence of anomalous heat conduction in carbon and 
boron-nitride nanotubes [5] (see also [6] for a review of recent experimental 
works). 

As shown by Matsuda and Ishii [7], in finite harmonic chains with uncor¬ 
related isotopic disorder most normal modes undergo Anderson localisation 
and practically do not contribute to heat transport. The situation is different 
for low-frequency modes, because their localisation length exceeds the size of 
the chain; these modes are essentially extended and carry energy across the 
system. In a random chain of N sites there are about y/N low-frequency de¬ 
localised modes [8]; their form depends crucially on the boundary conditions. 
Taking into account these facts, one can compute the heat flux with the help 
of the Matsuda-Ishii formula. After dividing the result by the temperature 
gradient, one eventually obtains that the conductivity scales as 

K oc AT", (1) 

with a = 1/2 for free boundary conditions and a = —1/2 for fixed boundary 
conditions [7, 9, 10]. The conductivity coefficient k,, therefore, depends on 
the size of the system: this implies that Fourier’s law cannot hold, because 
its validity rests on k being an intensive quantity. 

This phenomenon, known as anomalous heat transport, occurs in random 
harmonic chains with uncorrelated isotopic disorder. Significant changes of 
the behaviour of the thermal conductivity can be expected, however, when 
disorder exhibits spatial correlations. In fact, as made clear by the intensive 
research of the past fifteen years (see [11] and references therein), the localisa¬ 
tion properties of the eigenstates of low-dimensional systems can be strongly 
modified by imposing specific correlations to the disorder. Most studies of 
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ID models with correlated disorder have been focussed on the localisation of 
electronic states and electromagnetic or spin waves; less attention has been 
devoted to the effects of correlated disorder on the heat conduction in har¬ 
monic chains. Among the works centred on this topic, one should mention 
the analysis done in [12], where the authors numerically studied the diffusion 
of a localised energy input in random harmonic chains with long-range dis¬ 
order correlations. In Ref. [12] the sequence of random masses was obtained 
using the spectral method for the generation of a fractal Brownian motion. 
The same model was analysed in [13], where the structure of the normal 
modes and of the phonon spectrum was studied. 

Rather than focussing on long-range disorder correlations, in Ref. [14] the 
authors considered the effects of exponentially decaying correlations. They 
found that correlations of this kind lead to greater high-frequency phonon 
transmittance for short chains and to an increased thermal resistance for long 
chains. The impact of specihc long-range correlations was studied in Ref. [15] 
where it was shown how to delocalise a hnite fraction of the vibrational modes 
in high-frequency regions of the phononic spectrum. In this way one can 
obtain random chains with conductivity roughly proportional to the system 
size, i.e., k oc with a ~ 1, as occurs for harmonic chains without disorder. 

The previous works on random chains with correlated disorder did not 
establish the conditions required for the onset of normal heat conduction, 
but their results suggested that Fourier’s law might be recovered by imposing 
specihc correlations to the isotopic disorder. In this paper we show that this 
intuition is indeed correct and we present an approach which allows us to 
recover Fourier’s law in ID chains by controlling how k scales with the chain 
length N. 

Our idea goes as follows. Instead of using correlations to delocalise mid- 
to-high frequency modes, we now employ them to produce signihcant altera¬ 
tions of the fraction Nf,/N of low-frequency extended modes. In this way, by 
means of appropriate correlations, we can modify the scaling of the conduct¬ 
ivity with the system size. In particular, we show that, for hxed boundary 
conditions, specihc long-range correlations can produce a size-independent 
conductivity, thereby ensuring the validity of Fourier’s law even in ID har¬ 
monic chains. For free boundary conditions one does not obtain an intensive 
conductivity, but other long-range correlations can make the conductivity 
scale with the size of the chain as k ~ N^, with e positive but arbitrarily 
small. 
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2 The model 


We consider a ID harmonic chain of N atoms with nearest-neighbours inter¬ 
actions and the same elastic constant k thronghont the chain. The corres¬ 
ponding Hamiltonian has the form 




Here Pn and qn respectively represent the momentnm and the displacement 


from the equilibrinm position of the n-th atom of mass while D stands for 


the tridiagonal force-constant matrix with elements T)nm = ‘2k6n,m—k6n,m-i — 
kdn,m+i for n = 2,... , iV — 1 and m = 1, ..., iV. The hrst and the last rows of 
D determine the dynamical equation of the atoms at the extremities of the 
chain; their forms depend on the chosen boundary conditions. In the case of 
free boundary conditions the extreme atoms interact only with their single 
nearest neighbour, so that the hrst and last rows of the force matrix are 
Di,m = kSi^rn “ and = k6N,m “ k6N,m+i- For £xed boundary 

conditions, the outermost atoms of the chain are coupled with springs of 
elastic constant k not only to their nearest neighbours, but also to external 
walls of inhnite mass. Correspondingly, one has Di ^ = 2k6i^m — and 

~ ‘^kdpf^rn k6]\f^rn+l- 

Isotopic disorder enters the model (2) via random atomic masses 
which huctuate around a common average value 


M = (rrin). 


(3) 


In Eq. (3) and in what follows the symbol (■ • •) denotes the average over 
realisations of the isotopic disorder. We use the symbol 


(4) 


6mn = rUn — M 


for the huctuation of the n-th mass with respect to the average M. The 
strength of disorder is characterised by the variance of the huctuations 




(5) 


which we assume to be finite. Another key ingredient is the normalised binary 
correlator 



( 6 ) 
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whose specific form will be defined below according to the desidered proper¬ 
ties of the thermal conductivity. We stress that we do not need to assume 
the usual weak-disorder condition, 1. In the present work we are in¬ 

terested solely in the behaviour of the low-frequency modes, since only these 
modes come into play. For this reason Eqs. (5) and (6) represent a sufficient 
description of the statistical properties of the random masses. 

We remark that our model is characterised by isotopic disorder. One 
can consider also chains with random spring constants rather than random 
masses. However, the model with random masses is more appropriate for 
the study of thermal transport, since the isotopical doping is more feasible 
experimentally, see Ref. [6]. It should also be noted that random springs 
endow the model with off-diagonal disorder, whereas random masses generate 
purely diagonal disorder, which is mathematically easier to deal with. 

To study heat transfer, the extremities of the chain must be coupled to 
two heat baths of given temperatures T_ and T+. In the literature several 
models have been used for the baths. Two examples are: the Langevin 
representation, adopted in the Casher-Lebowitz model [9], and the semi- 
inhnite harmonic chain [16]. Recently, a generalised Langevin scheme was 
also considered in Ref. [4], demostrating the role of the spectral properties 
of the baths in shaping the transport properties of the chain. To make our 
approach more transparent, in this work we use the time-honoured Langevin 
representation of the baths. We would like to stress, however, that our results 
are equally valid for oscillators baths, as can be shown with the use of the 
formalism introduced by Dhar [4, 17]. 

In the Langevin scheme the action of the baths is mimicked by adding 
to the dynamical equations of the first and of the last atom an extra force 
composed of a dissipative term (proportional to the speed of the atom) and 
of a Gaussian white noise ^{t). The dynamical equations of the chain thus 
take the form 


Qn 

Pn 


Pn 

rrin 

72+1 


E D- 


7n=n—1 


JnN 


+ +-A 


I" +1 

Pn 

rriN 


- A 


mi 


(7) 


with n = 1,..., iV and go = <?w+i = 0. In Eq. (7), ^+(t) and +(t) represent 
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independent Gaussian white noises with zero average, ^±(t) = 0, and 


= 2XkBT+S{t-t') 
= 2XkBT_5 it - t') 
= 0 . 


( 8 ) 


We use the symbol (• ■ •) to represent the average over different realisations 
of the thermal noises. As can be seen from Eqs. (7) and (8), both baths are 
coupled with the same constant A to the chain. 

As is well known, in harmonic chains the heat flux is given by the sum 
of the modal fluxes. In the stationary regime and for weak coupling of the 
chain to the baths, the total flux is given by the Matsuda-Ishii formula [7] 


J = XkB (T+ 


T_ 


N 


r /-nI 

2 

//.a" 

1- 

1_ 


e)v 


i=i rriM 




1 2 


+ mi 


(*)' 

--N 


1 2 


(9) 


ii) 

where the symbols eh represents the n-th component of the i-th displacement 
eigenvector, as defined by the equation 

DeW=a;2g(q_ 

In Eq. (10), D is the rescaled force-const ant matrix with elements Dnm = 
Dnm/ ^/m^mm and the chain eigenfrequencies are labelled in ascending order: 

UJl < UJ2 < ... < Ci2jv 


3 Correlated disorder 


The Hamiltonian (2) defines the autonomous dynamics of the chain. Taking 
the Fourier transform of the dynamical equations with respect to time, one 
obtains 


/ cu A ^ 5mn 

<in+l + Qn-l + 4 I - I —— qn 

VcUmax/ M 
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Qn 


( 11 ) 


where 



( 12 ) 
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is the upper limit of the frequency spectrum for the homogeneous chain. In 
the absence of disorder, Eq. (11) reduces to the form 


(In+l + Q'n-l 
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qn 


( 13 ) 


which has solutions in the form of plane waves 


qn = C exp 


(14) 


with the wave number /i. Substituting the plane-wave solution in Eq. (13), 
one obtains the dispersion relation for the inhnite homogeneous chain 


UJ (/i) = Un 


1 — COS /i 


CUn 


. h 
sm — 
2 


( 15 ) 


The allowed values of the wavenumber p depend on the boundary conditions. 
One has 


with i 



Tii 

ivTT 

7r{i — 1) 
N 

In what 


for hxed boundary conditions 


( 16 ) 


for free boundary conditions 
follows we will use the short-hand notation 


Ui = 

When disorder is added, the normal modes suffer Anderson localisation, 
provided that the size of the chain is sufficiently large. We remark that the 
effective random potential felt by the Ath mode is given by 



f V 

yi^max / Af 


( 17 ) 


Taking into account Eqs. (4) and (5), it is easy to see that the effective 
potential (17) has zero average, = 0, and variance 


W‘T) = 16 


Ui 


OJr, 


M2' 


(18) 


We restrict our attention to the low-frequency modes, identihed by the con¬ 
dition 

CUj ^max- 
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For these modes the effective random potential is weak, as can be seen by 
considering that the deterministic part of the diagonal term in Eq. (11) is 
much larger than the random potential (17). In fact, when condition (19) is 
satisfied, one has 
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( 20 ) 


Note that condition (19) ensures that the relation (20) holds even if a/M ~ 1. 

For the low-frequency modes one can use the Hamiltonian map approach 
described in [15] and obtain the following expression for the inverse localisa¬ 
tion length 



2M2 


OOj 

^max 



( 21 ) 


In Eq. (21) the factor W represents the power spectrum of the random se¬ 
quence {6mn}, i.e., 

OO OO 

hF(p) = 1 -h 2 ^ x{l) cos (2/p) = ^ x(/) exp (i2p/). (22) 

1 = 1 l = — 00 


Following [7], we can now use formula (21) to estimate the number Ne of 
extended modes. In a finite chain, a mode can be considered extended if its 
localisation length lioc exceeds the chain length N, i.e., if the condition 


N < hoc(,(^i) = 


2M^ 




OJn 


/ hF (2()nj/Ci.)max) 


(23) 


is satisfied. Because low-frequency modes experience weak disorder, their 
eigenfrequencies stay close to the corresponding values for an homogeneous 
chain and one can still use Eq. (15). In the low-frequency region, the dis¬ 
persion law (15) can be linearised, so that in the end one can substitute in 
Eq. (23) the approximate expression 


Ui ~ Ua 


TTl 


( 24 ) 


with i N. Condition (23), supplemented by Eq. (24), leads to the con¬ 
clusion that the number of extended modes is obtained by solving the 
equation 




2M 


2N 




a 


W (ttNJN) ■ 


( 25 ) 













4 Emergence of normal conductivity 


To estimate the heat flux, we now evaluate the sum in Eq. (9) by taking into 
account that the main contribution is due to the extended modes, and 
that they have essentially the same form of those of an ordered chain. In this 
way, for iV 3> 1 one obtains 


Afe 


A/cb(T+— T_) . 2 

J ~ - . ^ > sm 


k=l 


MN 

XkB{T+-T_) /iVe 
M [ N 


Tik 


iV+ 1 


(26) 


for fixed boundary conditions and 

Ve 


J ~ 


XksiT^ — T„) 2 f 

- > cos —— 


MN 


k=l 


\kB{T+ - T_) iVe 


M 


N 


( 27 ) 


for free boundary conditions. Dividing the heat flux by the temperature 
gradient, one eventually arrives at the result that the conductivity scales as 


J iVg/iV^ for fixed boundary conditions 

iVe for free boundary conditions 


(28) 


For uncorrelated disorder W = 1 and with Eq. (25) we recover the well- 
known result [8] 

Ve ~ y/N. (29) 

Substituting this result in Eq. (28) one obtains that k I/a/ZV for fixed 
boundary conditions and that k ~ \/N for free boundary conditions. These 
remarkable results were established in the early 1970s by several authors [7, 
9, 16, 10] and led to the crucial conclusion that Fourier’s law cannot hold in 
ID disordered chains. Since then, this discovery has been the subject of an 
intense debate (see, e.g., [2] and references therein). 

The aim of this paper is to analyse whether spatial correlations of the 
disorder can produce an intensive conductivity and thus lead to normal heat 
conduction. We also want to explore the possibility to have other scaling laws 
for K as a function of N. To proceed in this direction one has to understand 
first how specific correlations of the isotopic disorder can change the number 
of low-frequency extended states and, therefore, alter the scaling behaviour 
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of the conductivity. Let us assume that in the long-wavelength limit the 
power spectrum ( 22 ) obeys a power law of the form 

hL(p) ~/i^ for /i —)■ 0. (30) 


Note that here f3 > —1; the values /? < —1 are forbidden because the power 
spectrum ( 22 ) must obey the normalisation condition 


t/2 


71 




(31) 


Substituting Eq. (30) into the condition (23), one obtains that the number 
of extended modes becomes 


Ne ~ 


(32) 


with 


1 + /5 

2 + 13 ' 


(33) 


Since /3 > — 1 , the values of 7 are restricted to the interval 0 < 7 < 1. 
After substituting the result (32) in Eq. (28), one can conclude that the 
conductivity obeys the power law ( 1 ) with scaling exponent a equal to 


a = 

a = 


/3-1 


13 + 2 
(3 + 3 


(3 + 2 


for hxed boundary conditions 
for free boundary conditions. 


(34) 


Eq. (34) represents the central result of this paper. It shows how disorder 
correlations can change the scaling of the conductivity with the size of the 
chain. 

Two special cases stand out for their interest. For hxed boundary condi¬ 
tions, one can consider the value (3 = 1 in Eq. (30). This implies that the 
power spectrum decreases linearly with the frequency in the limit a; —?■ 0 ; the 
corresponding binary correlator ( 6 ) decays with a power-law x(/) ~ 1 //^ for 
/ 3> 1. From Eq. (34), it is easy to see that for (3 = 1 the scaling exponent 
of the conductivity vanishes, so that the conductivity becomes independent 
of the size of the chain, k ~ iV°. The physical mechanism behind this res¬ 
ult is easy to understand. For (3 = 1, the inverse localisation length (21) 
scales as rather than being a quadratic function of the frequency 

as in the case of uncorrelated disorder. Therefore, the number of extended 
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modes is increased due to the chosen correlations, and one has Ar2/3 

instead of ~ The increase of summands in Eq. (26) compensates 

the 1/iV^ scaling of the displacement eigenvectors and produces an intensive 
conductivity. Therefore, Fourier’s law emerges in the asymptotic regime, i.e., 
for iV 3> 1. 

For free boundary conditions, a case of physical interest corresponds to 
a power spectrum diverging ioi cu — )■ 0 as W{u)) ~ This corresponds 

to /9 = — 1 + £, with £ being a small but positive parameter. In this case the 
correlation function decays extremely slowly over large distances: x(0 l/l^ 
for / 3> 1. The choice f3 = — 1 + e implies that the localisation length scales 
as l\oc ~ in the low-frequency region. Therefore the number of extended 
eigenmodes is now reduced with respect to the case of uncorrelated disorder; 
in fact, one has Ng ~ for e —)■ O’*". This reduction of the extended modes 
leads to the result 

n ~ ^ iVT (35) 

Since condition (31) requires e > 0, one cannot exactly obtain a normal 
conduction for free boundary conditions by manipulating the low-frequency 
behaviour of the power spectrum. However, for sufficiently small values of 
e, one can obtain a conductivity with a very weak dependence on the size of 
the chain. 

These theoretical predictions are conhrmed by the numerical data, as can 
be seen from Figs. 1 and 2, which show the scaling behaviour of the conduct¬ 
ivity for fixed and free boundary conditions. In numerical simulations we 
setM = l,/c = l,A = l, and = 0.2. We considered 100 different realisa¬ 
tions of the disorder and we computed the average value and the standard 
deviation of the conductivity over this ensemble. To generate random masses 
with the desired spatial correlations we have used the standard technique of 
filtering white-noise sequences [11]. 

For hxed boundary conditions, we chose a binary correlator of the form 


x (0 



^2 + (- 1 )'+^ 
1 - 16/2 


which corresponds to the power spectrum 


Wifi) 


niV2 + l) 
2^2 


sin 



(36) 


(37) 


The data represented in Fig. 1 clearly show how the conductivity flattens as 
N increases. 
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Figure 1: Conductivity k versus chain length N for hxed boundary condi¬ 
tions. Here we considered fluctuations of the random masses with the 
power spectrum (37) which, in the low-frequency region, behaves as Eq. (30) 
with (3 = 1. The parameters of the model used in simulations were: M = 1, 
k = 1, X = 1, and = 0.2. The points represent the mean value of the 
conductivity, obtained after averaging over 100 disorder realisations. Error 
bars represent the standard deviation of the numerical data with respect to 
the mean value. 


For free boundary conditions we generated a disorder with power spec¬ 
trum 

1—e 


fF(/i) = 


1 + (3 \7l 


3 

F) = 


-T 

2/i; 


(38) 


with e = 0.01 and e = 0.001. Fig. 2 demonstrates that, as e ^ 0+ and 
(3 = — 1 -f- £ approaches — I'*', the asymptotic behaviour of the conductivity 
becomes very close to the normal one. We remark that, as (3 moves closer 
to —1, the size of the ballistic region becomes larger and larger. This is to 
be expected, because the random masses become strongly correlated over 
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Figure 2: Conductivity k versus chain length N for free boundary conditions. 
The fluctuations of the random masses have the power spectrum (38) 
with (3 = 0.990 (corresponding to £ = 0.01) and (3 = 0.999 (corresponding to 
£ = 0.001). Here the parameters of the model were: M = 1, fc = 1, A = 1, and 
cr^ = 0.2. The points represent the mean value of the conductivity, obtained 
after averaging over 100 disorder realisations. Error bars are omitted because 
they are smaller than the symbols used to represent the mean values of k. 


very long distances; physically, this means that the chain is made up of large 
semi-homogeneous chunks. As N increases, the almost-normal behaviour is 
eventually reached. 

5 Conclusions 

We have analysed the heat conductivity k in harmonic chains with isotopic 
correlated disorder. As is known, in such chains normal heat conduction is 
not possible when the disorder is uncorrelated. Specihcally, the conductiv- 
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ity K, , instead of being independent of the length N of the chain, scales 
as K ~ A^“, where a depends on the boundary conditions. Our analytical 
approach shows that, with a proper choice of long-range correlations of the 
disorder, one can recover normal heat conduction for sufficiently large N. The 
length scale Nc for the onset of normal heat conduction, k iV° for N > iVc, 
strongly depends on the boundary conditions. Specihcally, for hxed bound¬ 
ary conditions one can obtain a conductivity strictly independent of N, while 
for free boundary conditions normal heat conduction can be achieved with an 
arbitrary precision. With our approach one can further show that, for trickier 
correlations, abnormal scaling laws can also emerge, such as k ~ A^/(lniV)‘^ 
with (5 = 1 for free boundary conditions and (5 = 3 for hxed boundary con¬ 
ditions. Although such a scaling occurs for extremely large values of N, this 
example shows that, in principle, one can get a conductivity which increases 
with the chain size regardless of the imposed boundary conditions. In gen¬ 
eral, our results demonstrate how the behaviour of the heat conductivity can 
be controlled by imposing specihc correlations in the disorder. This may 
prove important in future technological applications. 
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